###### RDF Code #######
install.packages(c("tidyr", "stringr", "plm", "foreign", "ggplot2", "haven", "ggplot2"))


library(tidyr)
library(stringr)
library(dplyr)
library(plm)
library(foreign)
library(ggplot2)
library(haven)
library(purrr)
library(sandwich)
library(lmtest)
library(readr)


setwd("Replication Files")


#Load RDF Data from https://www.pew.org/en/research-and-analysis/articles/2025/03/27/state-rainy-day-fund-growth-slowed-in-fiscal-2024
State_RDF <- read.csv("State_RDFs.csv")


#Load lcrps_merge
lcrps_merge <- readRDS("lcrps_merge")

#Subset State_RDF_data
State_RDF_sub <- State_RDF[,c(1:26)]

#Clean State_RDF data
State_RDF_wide <- pivot_longer(State_RDF_sub, cols=!State.Region, names_to="Year")
State_RDF_wide$Year <- str_replace(State_RDF_wide$Year,"Total.balances..days.","")
State_RDF_wide$value <- str_replace(State_RDF_wide$value,",","")
State_RDF_wide$value <- as.numeric(iconv(State_RDF_wide$value, 'utf-8', 'ascii', sub=''))
State_RDF_wide <-  State_RDF_wide %>% rename("State" = "State.Region")
State_RDF_wide$Year <- as.numeric(State_RDF_wide$Year)


#Merge Ranney with RDF
Ranney_RDF_merge <- State_RDF_wide %>% inner_join(lcrps_merge, by=c('State' = 'State', 'Year' = 'Year'))
Ranney_RDF_merge <- Ranney_RDF_merge %>% rename(days_gf = value)

saveRDS(Ranney_RDF_merge, "Ranney_RDF_merge")

#Tausanovitch and Warshaw Data
tw_data <- read.table('/Replication Files/aip_states_ideology_v2022a.tab',header=TRUE)

all_years <- 2008:2020

#Interpolate T-W data
interp_all <- tw_data %>%
  group_by(state) %>%
  summarise(
    Year = list(all_years),
    TW_Ideology = list(approx(presidential_year, mrp_ideology, xout = all_years)$y),
    .groups = "drop"
  ) %>%
  unnest(c(Year, TW_Ideology))


interp_all <- interp_all %>% rename(State = state)


#### PLOTS #####

#Folded Ranney 8 yr
ggplot(data = Ranney_RDF_merge %>% filter(days_gf < 400), mapping = aes(x = folded_ranney_8yrs, y = days_gf)) +
  geom_point()+
  geom_smooth(method="lm",se=F) +labs(colour= "Year", x ="Folded Ranney (8 yrs)", y = "Days of GF Expenses in RDF")+
  facet_wrap(~Year)+theme_bw()


##### REGRESSION CONTROLS ####

#New state econ covars
State_GDP <- read.csv("BEA_State_GDP.csv")
State_PI <- read.csv("BEA_Personal_Income.csv")
State_Unemp <- read.csv("State_Unemployment.csv")
Pop_00_09 <- read.csv("2000_2009_Pop.csv")
Pop_10_19 <- read.csv("2010_2019_Pop.csv")
Pop_20_24 <- read.csv("2020_2024_Pop.csv")


#Clean data
colnames(State_GDP) <- gsub("^X", "", colnames(State_GDP))
State_GDP <- State_GDP[,-1]

State_GDP <- State_GDP %>%
  rename(State = GeoName)

#Pivot state gdp to wide
State_GDP_Wide <- pivot_longer(State_GDP, cols=!State, names_to="Year")

colnames(State_GDP_Wide) <- c("State","Year","GDP")

#State PI
colnames(State_PI) <- gsub("^X", "", colnames(State_PI))
State_PI_Wide <- pivot_longer(State_PI, cols=!State, names_to="Year")

colnames(State_PI_Wide) <- c("State","Year","Personal Income")

#State Unemployment
#Average for each year
State_Unemp_Yearly <- State_Unemp %>%
  group_by(State, Year) %>%
  summarize(Avg.Unemp.Rate = mean(Rate, na.rm = TRUE)) %>%
  ungroup()


#Replace header
new_header <- as.character(Pop_00_09[1, ])
Pop_00_09 <- Pop_00_09[-1, ]
colnames(Pop_00_09) <- new_header

new_header <- as.character(Pop_10_19[1, ])
Pop_10_19 <- Pop_10_19[-1, ]
colnames(Pop_10_19) <- new_header

new_header <- as.character(Pop_20_24[1, ])
Pop_20_24 <- Pop_20_24[-1, ]
colnames(Pop_20_24) <- new_header

#Get rid of periods
Pop_00_09$State <- gsub("\\.", "", Pop_00_09$State)
Pop_10_19$State <- gsub("\\.", "", Pop_10_19$State)
Pop_20_24$State <- gsub("\\.", "", Pop_20_24$State)

#Merge all data frames
df_list <- list(Pop_00_09, Pop_10_19, Pop_20_24)

library(purrr)
State_Pop <- reduce(df_list, full_join, by = "State")

#Pivot data to long
State_Pop_Wide <- pivot_longer(State_Pop, cols=!State, names_to="Year")
colnames(State_Pop_Wide) <- c("State","Year","Population")



#Merge all the covariates data
df_list <- list(State_Pop_Wide, State_Unemp_Yearly, State_PI_Wide, State_GDP_Wide)
State_Pop_Wide$Year <- as.integer(State_Pop_Wide$Year)
State_Unemp_Yearly$Year <- as.integer(State_Unemp_Yearly$Year)
State_PI_Wide$Year <- as.integer(State_PI_Wide$Year)
State_GDP_Wide$Year <- as.integer(State_GDP_Wide$Year)

merge1 <- left_join(State_Unemp_Yearly,State_GDP_Wide, by = c("State","Year"))
merge2 <- left_join(merge1, State_PI_Wide, by = c("State","Year"))
State_Covs <- left_join(merge2, State_Pop_Wide, by = c("State","Year"))

State_Covs_2000s <- State_Covs %>% filter(2000 <= Year & Year <= 2024)

saveRDS(State_Covs_2000s,"State_Covs_2000s")

#State tax data
tax_rev <- read.csv('StateTaxRevenue.csv')
tax_rev_subset <- tax_rev[-1,c(1,3,4,5)]
names(tax_rev_subset) <- c("Year","Abbrev", "FY_Ending","Tot_Taxes")

tax_rev_subset$Abbrev <- substr(tax_rev_subset$Abbrev, 1, 2)
tax_rev_subset <- tax_rev_subset[-1,]

state_xwalk <- read.csv('StateName_Crosswalk.csv')
  
tax_rev_subset <- left_join(tax_rev_subset, state_xwalk, by = "Abbrev")

#Merge with Ranney data
Ranney_RDF_merge <- left_join(Ranney_RDF_merge, State_Covs_2000s, by = c("State", "Year")) 
Ranney_RDF_merge <- left_join(Ranney_RDF_merge, tax_rev_subset, by = c("State", "Year")) 

saveRDS(Ranney_RDF_merge, "Ranney_RDF_merge")
ranney_pdata <-  pdata.frame(Ranney_RDF_merge, index = c("State","Year"))

saveRDS(ranney_pdata,"ranney_pdata")

ranney_pdata <- readRDS('ranney_pdata')

#Read sle, merge with ranney_pdata, getting a rolling number of elections in the past 8 years, use that to calculate the mean party change
sle <- readRDS("sle")

sle$Abbrev <- toupper(sle$sab)
sle <- left_join(sle, state_xwalk, by = "Abbrev")
names(sle) <- c("sab", "Year","election_year", "Abbrev", "State", "FIPS")

Ranney_RDF_merge2 <- left_join(Ranney_RDF_merge, sle, by = c("State","Year"))

#Create a rolling count of elections
Ranney_RDF_merge2 <- Ranney_RDF_merge2 %>%
  arrange(State, Year) %>%
  group_by(State) %>%
  mutate(election_year = as.numeric(coalesce(election_year, 0)),
         elects_sum_roll = zoo::rollsum(election_year, k = 8, align = "right", fill = NA))

Ranney_RDF_merge2$mean_party_change_8yr <- Ranney_RDF_merge2$control_change_8yr/Ranney_RDF_merge2$elects_sum_roll

ranney_pdata2 <- pdata.frame(Ranney_RDF_merge2, index = c("State","Year"))

ranney_pdata2 <- ranney_pdata2[complete.cases(ranney_pdata2$mean_party_change_8yr),]


ranney_pdata <- ranney_pdata %>%
  mutate(
    GDP            = parse_number(as.character(GDP)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate)),
    Personal.Income        = parse_number(as.character(Personal.Income)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes))
  )

ranney_pdata2 <- ranney_pdata2 %>%
  mutate(
    GDP            = parse_number(as.character(GDP.x)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate.x)),
    Personal.Income        = parse_number(as.character(Personal.Income.x)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes)),
    mean_party_change_8yr      = parse_number(as.character(mean_party_change_8yr))
  )



#Merge ranney data with ideology data
Ranney_RDF_Ideology <- left_join(Ranney_RDF_merge, interp_all, by = c("State", "Year")) 
ranney_ideology_pdata <-  pdata.frame(Ranney_RDF_Ideology, index = c("State","Year"))

ranney_ideology_pdata <- ranney_ideology_pdata[complete.cases(ranney_ideology_pdata),]

ranney_ideology_pdata <- ranney_ideology_pdata %>%
  mutate(
    GDP            = parse_number(as.character(GDP.x)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate.x)),
    Personal.Income        = parse_number(as.character(Personal.Income.x)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes)),
    Population =  parse_number(as.character(Population.x))
  )

Ranney_RDF_Ideology2 <- left_join(Ranney_RDF_merge2, interp_all, by = c("State", "Year")) 

ranney_ideology_pdata2 <-  pdata.frame(Ranney_RDF_Ideology2, index = c("State","Year"))

ranney_ideology_pdata2 <- ranney_ideology_pdata2[complete.cases(ranney_ideology_pdata2),]

ranney_ideology_pdata2 <- ranney_ideology_pdata2 %>%
  mutate(
    GDP            = parse_number(as.character(GDP.x)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate.x)),
    Personal.Income        = parse_number(as.character(Personal.Income.x)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes)),
    Population =  parse_number(as.character(Population.x)),
    mean_party_change_8yr      = parse_number(as.character(mean_party_change_8yr))
  )

saveRDS(ranney_ideology_pdata,"ranney_ideology_pdata")

model2 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_4yrs + Avg.Unemp.Rate + GDP+ log(Tot_Taxes) +Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model2)


model3 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_6yrs +  Avg.Unemp.Rate + GDP+ log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model3)

model4 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_8yrs +  Avg.Unemp.Rate + GDP+log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model4)

model5 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_10yrs +  Avg.Unemp.Rate + GDP+log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model5)

model6 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_hl_index +  Avg.Unemp.Rate + GDP +log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "twosteps",
  transformation="d",
  collapse= TRUE,
  robust = TRUE
)

summary(model6)

model7 <- pgmm(
  days_gf ~ lag(days_gf,1) + mean_party_change_8yr +  Avg.Unemp.Rate + GDP+log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata2,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse= TRUE,
  robust = TRUE
)

summary(model7)

nobs(model7)

stargazer(model7)



library(stargazer)
stargazer(model2, model3, model4, model5, model6, model7,
          type = "latex",
          title = "Regression Results", 
          dep.var.labels = "RDF Size (Days of General Fund Support)", 
          covariate.labels = c("RDF Size (t-1)", "Folded Ranney (4 year)", "Folded Ranney (6 year)",
                               "Folded Ranney (8 year)","Folded Ranney (10 year)", "Hinchliffe and Lee Index", "Party Change",
                               "Avg. Unemployment Rate", "GDP", "Log Annual Tax Revenue", "Avg Personal Income"),
          ci = FALSE,  # You can't use ci=TRUE unless you're passing model objects that support it
          single.row = FALSE)

# Note: Stargazer is only reporting groups so the correct number of observations are
nobs(model2)
nobs(model3)
nobs(model4)
nobs(model5)
nobs(model6)
nobs(model7)


plm::mtest(model2, order =1)
plm::mtest(model3, order =1)
plm::mtest(model4, order =1)
plm::mtest(model5, order =1)
plm::mtest(model6, order =1)
plm::mtest(model7, order =1)


plm::mtest(model2, order =2)
plm::mtest(model3, order =2)
plm::mtest(model4, order =2)
plm::mtest(model5, order =2)
plm::mtest(model6, order =2)
plm::mtest(model7, order =2)

sargan(model2)
sargan(model3)
sargan(model4)
sargan(model5)
sargan(model6)
sargan(model7)


##### IDEOLOGY MODELS ######

ranney_ideology_pdata <- ranney_ideology_pdata[complete.cases(ranney_ideology_pdata),]
cor(ranney_ideology_pdata$ranney_8yrs, ranney_ideology_pdata$TW_Ideology)

model2I <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_4yrs + Avg.Unemp.Rate +GDP + log(Tot_Taxes)  + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model2I)

model3I <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_6yrs +Avg.Unemp.Rate +GDP + log(Tot_Taxes) + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model3I)

model4I <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_8yrs + Avg.Unemp.Rate +GDP + log(Tot_Taxes) + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model4I)

model5I <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_10yrs  + Avg.Unemp.Rate +GDP + log(Tot_Taxes) + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model5I)


model6I <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_hl_index  + Avg.Unemp.Rate +GDP + log(Tot_Taxes) + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model6I)

model7I <- pgmm(
  days_gf ~ lag(days_gf,1) + control_change_8yr  + log(Avg.Unemp.Rate) +GDP + log(Tot_Taxes) + TW_Ideology|
    lag(days_gf, 3:4),                  
  data = ranney_ideology_pdata2,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  robust= TRUE
)

summary(model7I)

stargazer(model2I, model3I, model4I, model5I,model6I,
          type = "latex",
          title = "Regression Results", 
          dep.var.labels = "RDF Size (Days of General Fund Support)", 
          covariate.labels = c("RDF Size (t-1)", "Folded Ranney (4 year)", "Folded Ranney (6 year)",
                               "Folded Ranney (8 year)","Folded Ranney (10 year)", "HL Index", "Party Change",
                               "Avg. Unemployment Rate", "GDP", "Log Annual Tax Revenue","Ideology"),
          omit.stat = c("LL", "ser", "f"),
          ci = FALSE,  # You can't use ci=TRUE unless you're passing model objects that support it
          single.row = FALSE)


# Number of observations
nobs(model2I)
nobs(model3I)
nobs(model4I)
nobs(model5I)
nobs(model6I)

#Tests
plm::mtest(model2I, order =1)
plm::mtest(model3I, order =1)
plm::mtest(model4I, order =1)
plm::mtest(model5I, order =1)
plm::mtest(model6I, order =1)


plm::mtest(model2I, order =2)
plm::mtest(model3I, order =2)
plm::mtest(model4I, order =2)
plm::mtest(model5I, order =2)
plm::mtest(model6I, order =2)

sargan(model2I)
sargan(model3I)
sargan(model4I)
sargan(model5I)
sargan(model6I)

mean(ranney_pdata$days_gf)


### WITHOUT WYOMING

Ranney_RDF_NoWY <- Ranney_RDF_merge %>% filter(State != "Wyoming")

ranney_pdata_NoWY <-  pdata.frame(Ranney_RDF_NoWY, index = c("State","Year"))


ranney_pdata_NoWY <- ranney_pdata_NoWY %>%
  mutate(
    GDP            = parse_number(as.character(GDP.x)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate.x)),
    Personal.Income        = parse_number(as.character(Personal.Income.x)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes)),
    Population =  parse_number(as.character(Population.x))
  )

Ranney_RDF_NoWY2 <- Ranney_RDF_merge2 %>% filter(State != "Wyoming")
Ranney_RDF_NoWY2 <- Ranney_RDF_NoWY2[complete.cases(Ranney_RDF_NoWY2),]
ranney_pdata2 <-  pdata.frame(Ranney_RDF_NoWY2, index = c("State","Year"))

ranney_pdata2 <- ranney_pdata2 %>%
  mutate(
    GDP            = parse_number(as.character(GDP.x)),
    Avg.Unemp.Rate = parse_number(as.character(Avg.Unemp.Rate.x)),
    Personal.Income        = parse_number(as.character(Personal.Income.x)),
    Tot_Taxes      = parse_number(as.character(Tot_Taxes)),
    Population =  parse_number(as.character(Population.x)),
    mean_party_change_8yr = parse_number(as.character(mean_party_change_8yr))
  )


model2 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_4yrs + Avg.Unemp.Rate + GDP+ log(Tot_Taxes) +Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata_NoWY,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model2)

sargan(model2,"twostep")



model3 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_6yrs +  Avg.Unemp.Rate + GDP+ log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata_NoWY,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model3)

model4 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_8yrs +  Avg.Unemp.Rate + GDP+log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata_NoWY,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model4)

model5 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_ranney_10yrs +  Avg.Unemp.Rate + GDP+log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata_NoWY,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse = TRUE,
  robust = TRUE
)

summary(model5)

model6 <- pgmm(
  days_gf ~ lag(days_gf,1) + folded_hl_index +  Avg.Unemp.Rate + GDP +log(Tot_Taxes)+Personal.Income|
    lag(days_gf, 3:4),                  
  data = ranney_pdata_NoWY,
  index = c("State","Year"),
  effect = "twoways",
  model = "onestep",
  transformation = "d",
  collapse= TRUE,
  robust = TRUE
)

summary(model6)


stargazer(model2, model3, model4, model5, model6,
          type = "latex",
          title = "Regression Results", 
          dep.var.labels = "RDF Size (Days of General Fund Support)", 
          covariate.labels = c("RDF Size (t-1)", "Folded Ranney (4 year)", "Folded Ranney (6 year)",
                               "Folded Ranney (8 year)","Folded Ranney (10 year)", "Hinchliffe and Lee Index", "Party Change",
                               "Avg. Unemployment Rate", "GDP", "Log Annual Tax Revenue", "Avg Personal Income"),
          ci = FALSE,  # You can't use ci=TRUE unless you're passing model objects that support it
          single.row = FALSE)

nobs(model2)
nobs(model3)
nobs(model4)
nobs(model5)
nobs(model6)



plm::mtest(model2, order =1)
plm::mtest(model3, order =1)
plm::mtest(model4, order =1)
plm::mtest(model5, order =1)
plm::mtest(model6, order =1)


plm::mtest(model2, order =2)
plm::mtest(model3, order =2)
plm::mtest(model4, order =2)
plm::mtest(model5, order =2)
plm::mtest(model6, order =2)

sargan(model2)
sargan(model3)
sargan(model4)
sargan(model5)
sargan(model6)



